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Abstract. We compare the luminosity function and rate inferred from the BATSE short hard bursts (SHBs) 
peak flux distribution with the redshift and luminosity distributions of SHBs observed by Sm/i/HETE II. While 

■ the Suii/t/HETE II SHB sample is incompatible with SHB population that follows the star formation rate, it is 
_^ ' compatible with a SHB rate that reflects a distribution of delay times after the SFR. This would be the case if SHBs 

are associated with binary neutron star mergers. The available data allows, however, different interpretations. For 
example, a population whose rate is independent of the redshift fits the data very well. The implied SHB rates 
that we find range from ~ 8 to ~ 30/i7 Gpc -3 yr _1 . This rate, which is comparable to the rate of neutron star 
mergers estimated from statistics of binary pulsars, is a much higher rate than what was previously estimated. 
We stress that our analysis, which is based on observed short hard burst is limited to bursts with luminosities 
above 10 49 erg/sec. Weaker burst may exist but if so they are hardly detected by BATSE of Swift and hence their 
rate is very weakly constrained by current observations. 

f--^ . Key words. cosmology:observations-gamma rays:bursts-gravitational radiation 

^ : 

1. Introduction 

2 ' It has been known since the early nineties that gamma-ray bursts (GRBs) are divided to two subgroups of long and 
short according to their duration: Tgo<2sec (Kouveliotou et al. 1993). Short bursts are also harder than long ones 
CS ■ (Dezalay et al. 1996, Kouveliotou et al. 1996, Qin et al. 2000) and hence are denoted Short Hard Bursts (SHBs). Our 
^ understanding of long burst and their association with stellar collapse followed from the discovery in 1997 of GRB 
^ ■ afterglow and the subsequent identification of host galaxies, redshift measurements and even detection of associated 
Supernovae. However, until recently no afterglow was detected from any short burst and those remained as mysterious 
as ever. 

This situation has changed with the detection of X-ray afterglow from several short bursts by Swift (Gehrels et al. 
2005, Romano et al. 2005) and by HETE II (Villasenor et al. 2005). In some cases optical (Covino et al. 2005, Fox et al. 
2005, Bloom et al., 2005, Price et al. 2005,Jensen et al. 2005, Hjorth et al. 2005, Gladders et al. 2005, Castro-Tiraldo 
et al. 2005, Gal- Yam et al. 2005, Cobb et al. 2005, Wiersema et al. 2005) and radio (Cameron & Frail 2005, Berger 
2005) afterglow was detected as well. This has lead to identification of host galaxies and to redshift measurements. 
While the current sample is very small several features emerge. First, unlike long GRBs that take place in galaxies 
with young stellar population SHBs take place also in elliptical galaxies in which the stellar population is older. In 
this they behave like type la Supernovae. The redshift and peak (isotropic equivalent) luminosity distributions of the 
five short bursts (see Table ^) indicate that the observed SHB population is significantly nearer than the observed 
long burst population. This feature was expected (Piran 1994; Katz & Canel 1996; Tavani 1998, Guetta & Piran 2005 
denoted hereafter GP05) as the (V/V max ) = 0.39 ± 0.02 of the BATSE short burst population was significantly larger 
(and closer to the Euclidian value of 0.5) than the one of long bursts ((V/V max ) = 0.29 ± 0.01, Guetta, Piran & 
Waxman 2004 denoted hereafter GPW). 

Recently GP05 have estimated the luminosity function and formation rate of SHBs from the BATSE peak flux 
distribution. These two quantities are fundamental to understanding the nature of these objects. The observed flux 
distribution is a convolution of these two unknown functions, so it is impossible to determine both functions without 
additional information. Already in 1995 Cohen and Piran (1995) have shown that the observed BATSE flux distribution 
can be fitted with very different luminosity functions depending on the choice of the GRB rate. GP05 have shown 
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0.16 
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L 7 ,i 3O /10 bl erg/sec 


0.14 
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0.17 


1.9 
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Table 1. The Swift/KETE II current sample of SHBs with a known rcdshift 



that the distribution is compatible with cither a population of sources that follow the SFR (like long bursts) or with 
a population that lags after the SFR (Piran, 1992, Ando 2004). There are several reasons to expect (see e.g. Narayan, 
Piran & Kumar 2001) that SHBs may be linked to binary neutron star mergers (Eichlcr ct al., 1989). In such a case the 
SHB rate is given by the convolution of the star formation rate with the distribution P m (r) of the merging time delays 
t of the binary system. These delays reflect the time it takes to the system to merge due to emission of gravitational 
radiation. 

As BATSE is less sensitive to short bursts than to long ones (Mao, Narayan & Piran, 1994), even an intrinsic 
SHB distribution that follow the SFR gives rise to an observed distribution that is nearer to us. Still a delayed SFR 
distribution (that is intrinsically nearer) gives rise to even nearer observed distribution (GP05). Therefore the recent 
observed rcdshift distribution of SHBs favors the delayed model and hence the merger scenario. Still the question 
was posed whether the predicted (GP05) observed distribution is consistent with the current sample. Gal Yam et al., 
(2005) suggested that the distributions are inconsistent and hence the suggested delayed SFR model is ruled out. We 
re-examine the situation here and we show that while a delayed distribution with "maximal likelihood" parameters are 
indeed ruled out, a delayed distribution with parameters within la from the best fit parameters cannot be ruled out 
with the current SHB rcdshift distribution. We discuss the implications of this result to GRBs and to binary Neutron 
star mergers and to the detection of gravitational radiation from such mergers. 



2. The luminosity function of the BATSE SHB sample 

Our data set and methodology follow GP05. We consider all the SHBs detected while the BATSE onboard trigger 
(Paciesas et al. 1999) was set for 5.5a over background in at least two detectors in the energy range 50-300keV. These 
constitute a group of 194 bursts. We assume, following a physical model, a rate of bursts. We then search for a best 
fit luminosity function. Using this luminosity function we calculate the expected distribution of observed redshifts and 
we compare it with the present data. 

We consider the following cosmological rates: 

— (i) A rate that follows the SFR. We do not expect that this reflects the rate of SHBs but we include this case for 
comparison. 

— (ii) A rate that follows the NS-NS merger rate. This rate depends on the formation rate of NS binaries, that 
one can safely assume follow that SFR, and on the distribution of merging time delays. This, in turn, depends 
on the distribution of initial orbital separation a between the two stars (r oc a 4 ) and on the distribution of 
initial eccentricities. Both are unknown. From the coalescence time distribution of six double neutron star binaries 
(Champion et al. 2004) it seems that P(log(r))dlog(r) ~ const, implying P m (r) oc 1/t, in agreement with the 
suggestion by Piran (1992). Therefore our best guess scenario is a SBH rate that follows the SFR with a logarithmic 
time delay distribution. In this case the normalization of P(t) is such that P(t) ^ only for 20Myr < r < age of 
the universe. Obviously delays longer than the age of the universe do not add events. 

— (iii) A rate that follows the SFR with a delay distribution P(r)dr ~ const. In this case the normalization of P(t) 
is such that P(t) ^ only for < r < age of the universe. 

— (iv) A constant rate (which is independent of rcdshift.). 

For the SFR needed in distributions (i— iii) we employ the SF2 model of Porciani & Madau (2001): 

_ 23exp(3.4z)[r»M(l + z) 3 + ^(l + z) 2 + fiA] 1/2 m 
RSF2{Z) P0 (exp(3.4z) + 22)(l + z)3/2 ' (1) 

where po is the present GRB rate and SlM,A,k are the present cosmological parameters. In models (ii) and (iii) the rate 
of SHBs is given by: 

,t(z) 

Rshb{z) = C 1 I dTR SF2 (t-T)P m (T), (2) 
Jo 

where P m {r) is the distribution of the merging time delays r and C\ a normalization constant. 
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Table 2. Best fit parameters Rate(z = 0) , L* , a and /3 and their la confidence levels for models (i)-(iv). Also 
shown are the goodness of the fit and the KS probability that the five bursts with a known redshift arise from this 
distribution. We show two results for KS tests one with GRB 0508132 at z = 0.7 and the other at z = 1.8. Case ii CT 
corresponds to case ii with an L* value lower by la than the best fit one. Other parameters have been best fitted for 
this fixed number. 

In Fig. 1 Cmax is the count rate in the second brightest illuminated detector and C m in is the minimum detectable rate. 



Following Schmidt (2001), GPW and GP05 we consider a broken power law peak luminosity function with lower 
and upper limits, 1/Ai and A 2 , respectively: 

^ , rU1 r „ r f (L/L*)- a L*/Ai < L < L* ... 
Z o (L)dlogL = C dlogL^ {L ' /L J ) _ f3 L J <L<AiL . , (3) 

where Co is a normalization constant. This is the "isotropic-equivalent" luminosity function, i.e. it does not include a 
correction factor due to the fact that GRBs are beamed. In our calculations we approximate, following Schmidt (2001) 
the typical effective spectral index in the observed range of 20 or 50kcV to 300 keV as —1.1 (N(E) oc E~ 1A ). 

Following GP05 we use Ai i2 = (30, 100). Both values are chosen in such a way that even if there are bursts less 
luminous than L* / Ai or more luminous than A 2 L* they will be only very few (less than about 1%) of the observed 
bursts outside the range (L * /Ai, L*A 2 ). Bursts above L*A 2 are very bright and are detected to very large distances. 
However, there are very few bursts above L*A 2 . IncreasingA2 does not add a significant number of bursts (observed 
or not) and this does not change the results. In particular it does not change the overall rate. Ai is more subtle. The 
luminosity function increases rapidly with decreasing luminosity. Thus, a decrease in Ai will have a strong effect on the 
overall rate of short GRBs. However, it is important to realize that most of the bursts below L*/A\ are undetectable by 
current detectors, unless they are extremely nearby. Even if the luminosity function continues all the way to zero. This 
will increase enormously the over all rate of the bursts (which will in fact diverge in this extreme example) however 
most of these additional weak bursts will be undetected and the total number of detected bursts won't increase. For 
example, in our ii CT model only 1% of the BATSE detected bursts would have, in this extreme case, a luminosity below 
L * /A\. Clearly, the BATSE data does not constrain this part of the luminosity phase space and the question what 
is the number of such weak bursts is open, at least as far as the observed flux distribution is concerned. 

Comparing the predicted distribution with the one observed by BATSE we obtain, using a minimum \ 2 method, 
the best fit parameters of each model and their standard deviation (see Table 2). In order to assign a x 2 value we 
divided the 193 bursts in the BATSE catalog into 20 bins of equal size according to their value of P, the peak photon 
flux. Since the overall normalization is an additional free parameter the number of degree of freedom is 20-4=16. The 
results are presented in Table 2 and in Figs. 1 and 2 which depict comparisons between the predicted integrated and 
differential distributions and the observed one (Using the best fit parameters for the luminosity function.) For each 
parameter the la range is marginalized over the variation in all other parameters. The case iio- corresponds to the 
case ii with an L* value lower by la than the best fit one and we show also this model in Fig. 1. Below this value of 
L* the quality of the fit decreases very quickly. The possibility to decrease L* by a factor 10 and still get a reasonable 
fit is a crucial point of our analysis. In fact this is the key to making the BATSE results compatible with the low z 
distribution observed by Swift as we show in the next section. This is also the reason why the current estimates give 
much higher rate. 



3. A Comparison with the current Swift-HETE II SHB sample 

We can derive now the expected redshift distribution of the observed bursts' population in the different model.: 
N(z) = W*)^*) $ o( L)dlogL , (4) 



1 + z dz 



where £ max = A 2 L* = lOOi* and £ m in(-Piim) is the luminosity corresponding to minimum peak flux Pu m for a burst 
at redshift z. This minimal peak flux corresponds to the detector's threshold. For BATSE Pu m ~ 1 ph/cm 2 /sec for 
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Fig. 1. The predicted logN-log(P/Pii m ) distribution vs. the observed logN-log(C max /C m i n ) taken from the BATSE 
catalog for the best fit values of a,/3 and L* with models (i)-(iv). Also shown is the curve where the L* value is lower 
by la than the best fit one for the case (ii) we call this model (ii^) 



short bursts. The triggering algorithm for Swift is rather complicated but it can be approximated by a comparable 
minimal rate: Pu m (Swift) ~ 1 ph/cm 2 /sec. 

Fig-Eldepicts the SHB rate density as a function of redshift for the different r distributions and compare the results 
with a distribution that follows the SFR (SF2). As expected the time delay increases the number of short bursts at 
low redshift. In particular there is a dramatic increase if P(r)dr ~ const. Figs. 0] and [5] depict the expected observed 
(differential and integrated respectively) redshift distributions of SHBs in the different models. 

As expected, a distribution that follows the SFR, case (i), is ruled out by a KS test with the current five bursts. 
A KS test indicates that such a distribution is ruled out (the p- value is smaller than 1%). This is not surprising as 
other indications, such as the association of some SHBs with elliptical galaxies suggest that SHBs are not associated 
with young stellar populations. 

A distribution that follows the SFR with a constant logaritmic delay distribution, case (ii), is more interesting. 
A priori this was our favorite distribution as a logarithmic distribution of separations between the two members of a 
NS-binary (or a BH-NS binary) would lead to a logarithmic distribution of time delays (Piran 1992). At first glance 
this distribution is marginally consistent with the data. A KS test suggests that the probability that the observed data 
arises from this distribution is ~ 5-6% (see also Gal- Yam., 2005). The observed bursts are nearer (lower redshift) than 
expected from this distribution. It is interesting that the addition of a fifth burst (GRB 051221A) that was detected 
after a first version of this paper was submitted improved slightly the KS test value of this model. Additionally if we 
use the Rowan- Robinson SFR (see GP05), rather than SF2 of Porciani-Madau, the KS test for this model is 10% (for 
cither z=0.7 or z-1.8). 

While it is possible that the real distribution of time delays is not logarithmically constant there are several other 
possible explanations for this result. First we turn to the data and realize that selection effects that determine which 
SHB is detected and localized with sufficient accuracy to allow redshift determination are not clear. It is possible that 
we are dealing with small number statistics. It is also possible that the afterglows of these five localized bursts are 
brighter than the afterglow of a typical more distant bursts and this has influenced the sample. A second possibility 
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Fig. 2. The predicted differential distribution, n(P/Plim), vs. the observed n(C max /C m i n ) taken from the BATSE 
catalog for the best fit values of a, ft and L* with models (i)-(iv) and (iio-). 

is that the current data is a good sample of the SHB distribution but the "best fit" parameters estimated using the 
BATSE SHB population are, due to a statistical fluctuation, slightly offset. For example, we have considered in case 
(iicr) a distribution whose typical luminosity, L* , is one a away from the maximal likelihood value. This distribution 
that is consistent with the BATSE SHB sample is not ruled out by the current sample of SHBs with a known redshift. 
The p value of the KS test is 0.1. 

As an example for the flexibility of the data we have considered two other time delay distributions. Case (iii) in 
which the time delay distribution is uniform and case (iv) in which the overall SHB rate is constant in z. Both cases 
are compatible with the BATSE SHB distribution and with the sample of SHBs with a known redshift (The KS p 
values are 0.8 and 0.4 respectively.). This result is not surprising. The BATSE peak flux distribution depends on two 
unknown functions, the rate and the luminosity function. There is enough freedom to chose one function (the rate) 
and fit for the other. The sample of SHBs with known redshifts peaks at a rather low redshift and both distributions 
considered above push the intrinsic SHB rate to lower redshifts (as compared to the SFR). 

In all cases that are compatible with the five bursts with a known redshift, the intrinsic SHB rate is pushed 
towards lower redshifts and the present rates (needed to fit the observed BATSE rate) are larger by almost two orders 
of magnitude than those estimates earlier assuming that SBHs follow the SFR or that they follow it with a logarithmic 
delay (with the best fit parameters). These rate are ~ 30, ~ 8 and ~ 10/i7 Gpc~ 3 yr~ 1 for cases (iii), (iv) and (ii CT ) 
respectively. The corresponding "typical" luminosities, L*, ranges from 0.1 to 0.7 x 10 51 erg/sec. 

We have also considered, but did not present here for lack of space delay distributions that are peaked at a 
given delay (such as a log-normal distribution). A narrow distribution results simply in an SFR with a fixed delay. 
If the typical delay is long enough it will push the intrinsic SHB rate to low enough redshift. However, the current 
distribution of merging time of binary pulsars (Champion et al., 2004, GP05) indicates that the distribution is wide. 
A wide distribution, for example, a distribution with a typical delay of 10 9 years and a spread of two decades results 
in an overall rate that is very similar to one resulting from a logarithmic distribution. 

Before turning to the implications of these results we should, however, re-examine the fit of the data sets to the 
model. One way to do so is to compare the expected 2 dimensional distribution of bursts in the (z, log(L)) plane with 
the observed Swift/HETE II bursts. Fig El depicts such a comparison for our best model (iii). Since GRB 050709 was 
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Fig. 3. The rate of SHBs for the different models (i)-(iv). Note that the event rate of models ii and ii CT are the same 
as the difference between the two models is in the luminosity function and not in the event rate. 



determined by HETE we have reproduced the contours considering both the Swift and HETE sensitivity. While two of 
the bursts are just where we expect them to be the others are in low probability regions. Similar behavior arises for the 
other models as well. One can apply a sophisticated 2D KS test (Bahcall et al., 1987) to estimate the probability that 
the data and the models are consistent but given the small number of bursts we feel that this is not necessary at this 
stage. It is sufficient to say that an eye inspection suggests an inconsistency between all models and the Swift/HETE 
II data. 

There are several ways in which one can explain this apparent inconsistency. Once more we can turn to small 
number statistics and argue that it is difficult to make far reaching conclusions from such a few bursts. However, it 
is also possible that the data indicates a real discrepancy between the BATSE data (on which the models are based) 
and the Swift /HETE II data (to which the models are compared). This is not surprising as the triggering algorithm 
as well as the energy bands in which different bursts are observed in these two experiments are drastically different. It 
is possible for example that we are underestimating the limiting peak flux needed for detection by Swift. An increase 
in the estimated threshold of Swift would have improved the fit. Another issue that could influence the comparison 
is the existence of selection effects that influence the delegability of redshift of SHBs. Such selection effects could, of 
course, influence the sample of observed bursts with known redshift. 

It is also possible that these results indicate that there is something new in the data. A glance at Fig. 2 (the 
differential distribution of observed fluxes of short BATSE bursts) indicates a deep at an intermediate- high peak flux 
level, at P = 6photons/cm 2 /sec - about 6 times the minimal detectable flux. This is a flux level in which there is no 
observational reason why BATSE should have missed bursts. Still a marginally significant deep is there. It is possible 
that this deep reflects a real phenomenon and that there are two populations of SHBs? One that gives rise to the 



Guetta & Piran: Short duration GRBs 



7 




Z 



Fig. 4. The expected observed differential redshift distributions of SHBs for models (i)-(iv) and (iio-). The known 
redshifts of five SHBs are shown as vertical lines. The two small vertical lines corresponds to z=0.7 and z=1.8 



highest fluxes (at a level of ~ 10photons/cm 2 /sec) and another that gives rise to the low flux ones. A related conclusion 
has been reached recently using independent reasoning by Tanvir et al., (2005). 

4. Conclusions and Implications 

We have repeated the analysis of fitting the BATSE SHB data to a model of the luminosities and rates distributions. 
We stress that the fit to the BATSE data is similar to the one presented in GP05. Our best fit logarithmic distribution 
model is similar to the best fit logarithmic model presented in GP05. Our results and in particular the estimated rates 
seem very different from those obtained in GP05 but a glance in table 2 shows that the la spread in the rate is indeed 
very large! A main new ingredient of this work is the fact that we consider several other models. All fit the BATSE 
data equally well. We confirm our earlier finding that the BATSE data allows a lot of flexibility in the combination 
of the rates and luminosities. In GP05 we have shown that the BATSE data is compatible with a distribution that 
follows the SFR and one that has a logarithmic time lag. Here we consider additional distributions as well. 

A second new ingredient of this work is the comparison of the best fit models to the small sample of five Swift/HETE 
II SHBs. The Swift/HETE II data gives a new constraint. This constraint favors a population of SHBs with a lower 
intrinsic luminosity and hence a nearer observed redshift distribution. If this interpretation is correct it implies a 
significantly higher local SHB rate - a factor 50 or so higher than earlier estimates. We stress again that this new 
result was within the la error of the model presented in GP05. The difference between the conclusions in GP05 and 
the conclusions of this paper arises because of the new observations of Swift that show that the SHBs are nearer than 
what was expected before. This can arise due to one of the following possibilities: 1) The logarithmic time delay model 
is wrong. 2) The model is right but the best fit parameters are la away from the real ones or that the small Swift 
sample is misleading. 3) The Swift SHB set have significantly different selection effects than the BATSE SHB set. 



8 



Guetta & Piran: Short duration GRBs 



Empirical CDF 



1 

0.9 
0.8 
0.7 
0.6 
So.5 
0.4 
0.3 
0.2 
0.1 




P=1/t (mod ii) 
P=1/t 1-o (mod ii) 
P=K(mod iii) 
no z ev. (mod iv) 
SF2-sfr (mod i) 
obs. (z=1.8) 
obs. (z=0.7) 







0.5 



1.5 



2.5 



3.5 



Fig. 5. A comparison between the expected integrated observed redshift distributions of SHBs for models (i)-(iv) and 
(ii CT ) and the distribution of known redshifts of SHBs. 



Such a situation would arise if there are two populations of short bursts that are detected at different combinations 
by Swift and by BATSE. 

Provided that the basic model is correct and we are not mislead by statistical (small numbers), observational 
(selection effects and threshold estimates) of intrinsic (two SHB population) factors we can proceed and compare the 
inferred SHB rate with the observationally inferred rate of NS-NS mergers in our galaxy (Phinney, 1992, Narayan, 
Piran & Shemi 1992). This rate was recently reevaluated with the discovery of PSR J1829+2456 to be rather large as 
80^gg°/Myr. Although the estimate contains a fair amount of uncertainty (Kalogera et al. 2004). If we assume that this 
rate is typical and that the number density of galaxies is ~ 10~ 2 /Mpc 3 , we find a merger rate of 800±§gg°/Gpc 3 /yr. 
Recently Berger et al. (2005) have derived a beaming factor of 30-50 for short bursts. This rate implies a total merger 
rate of ~ 240 — 1500/Gpc 3 /yr for the three cases (iii), (iv) and (iio-). The agreement between the completely different 
estimates is surprising and could be completely coincidental as both estimates are based on very few events. 

If correct these estimates are excellent news for gravitational radiation searches, for which neutron star mergers 
are prime targets. They imply that the recently updates high merger rate, that depends mostly on one object, PSR 
J0737-3039, is valid. These estimate implies one merger event within ~ 70Mpc per year and one merger accompanied 
with a SHB within ~ 230Mpc. These ranges are almost within the capability of LIGO I an certainly within the 
capability of LIGO II. In these estimates we consider a luminosity function with a sharp cut off at a few xl0 49 erg/sec. 
Bursts with lower peak luminosity are practically undetectable by current detectors. Such bursts, if exist, constitute a 
very small fraction (a few percent at most) of the observed burst population. Therefore, using the current GRB data 
we cannot rule our or verify their existence. From this point of view one should consider our estimates of the rates of 
SHBs and hence the rate of neutron star mergers as lower limits. It is possible that the actual rate is much higher. If 
true this could be tested by LIGO I and II within the next few years. 

This research was supported by the US-Israel BSF and by the Schwarzmann university chair (TP). 
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Fig. 6. The two dimensional probability distribution of expected redshift and luminosity for model (iii) considering 
the sensitivity of Swift (black) and HETE (green). Contour lines are 0.9,0.8... 0.1 of the maximum. Also marked are 
the five SHBs with a known redshift. 
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